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Abstract 

We derive new formulas for the fundamental solutions of slow, viscous flow, governed by the 
Stokes equations, in a half-space. They are simpler than the classical representations obtained by 
Blake and collaborators, and can be efficiently implemented using existing fast solvers libraries. 
We show, for example, that the velocity field induced by a Stokeslet can be annihilated on the 
boundary (to establish a zero slip condition) using a single reflected Stokeslet combined with 
a single Papkovich-Neuber potential that involves only a scalar harmonic function. The new 
representation has a physically intuitive interpretation. 


1 Introduction 

Viscous flow of passive and active suspensions in the presence of an infinite planar boundary 
is an important physical model in many areas of science and engineering. It serves as a useful 
paradigm for understanding the effect of confined geometries on the macroscopic flow behavior 
of particulate flows, for example, that of bacterial propulsion, cellular blood flow and colloidal 
suspensions, |15l [TJ [3l [23l | 2 T] . In problems where the Reynolds number is low, the ambient fluid is 
governed by the Stokes equations: 

^ Au(x) = Vp(x), V • u(x) = 0, (1) 

where /r is the fluid viscosity, u(x) = (ui(x), M 2 (x), U 3 (x)) is the velocity of the fluid, and p(x) is 
the pressure. Assuming the plane wall is located at X 3 = 0 and that the flow velocity decays in the 
far held, 

u(x) —)■ 0 , as |x| —)■ 00 , 
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the no-slip boundary condition is 


u 


X3=0 


= 0 . 


( 2 ) 


Boundary integral methods are particularly well-suited for problems of this kind since they 
discretize the domain boundary alone (resulting in many fewer degrees of freedom) and impose the 
decay and the no-slip conditions exactly. Moreover, they avoid the need for artificial truncation of 
the computational domain and can be solved rapidly and with great accuracy using fast algorithms 
such as the fast multipole method (FMM) and high order accurate quadrature rules. 

In order to reformulate the Stokes equations as a boundary integral equation, however, one 
needs to have access to the Green’s functions for the half-space [mm]. A now classical approach 
to constructing this Green’s function is due to Blake and others mia min [21 [201 [ 29 ]. Unfortunately, 
the resulting formulas are rather complicated, making them somewhat difficult to implement. 

Here, we show that a much simpler alternative to the Blake solution can be obtained by combin¬ 
ing a free space image, which annihilates the tangential components of velocity, with a Papkovich- 
Neuber potential mm which annihilates the normal component. 

The present paper is organized as follows. We discuss Papkovich-Neuber potentials and the 
standard fundamental solutions for the Stokes equations in section (See, for example, [HI [HI 
[HlEO]). We also review Blake’s formula for the Stokeslet in a half-space. In section we derive 
the image structures for Stokeslets, stresslets, rotlets, and Stokes doublets. 

In the appendices, we provide the analogous formulas for two-dimensional half-space Stokes 
kernels and for some problems of linear elasticity. 


2 Fundamental solutions, the Papkovich-Neuber representation, 
and Blake’s formulas 


Before discussing the various standard fundamental solutions for the Stokes equations in free-space, 
we introduce the Papkovich-Neuber representation, originally developed in |19[ [T8] for problems of 
linear elasticity. Without loss of generality, we assume that the fluid viscosity /r = 1 in the rest of 
the paper. 

Definition 2.1. Let 4>{x),x G IR^ be a harmonic function. Then the induced Papkovich-Neuber 
representation is defined to be the paired vector field u and scalar field p given by: 


u(a;) = X 3 V<f>{x) 


0 

0 

(fix) 


p{x) = 2 


d(f{x) 

dx3 


( 3 ) 


It is straightforward to verify that (u,p) satisfy the Stokes equations Q (with ^ = 1). 

Suppose now that a force vector f = (fi, f 2 , f 3 ) is applied to a viscous fluid at a point y. Then, 
it is well-known that the induced velocity and pressure can be computed using the Stokeslet (the 
single layer kernel): 


S'ii(x,y) 


1 

Svr 


Sjj {xj - yi){xj - Vj) 

x-y| |x-yP 


LJ = 1,2,3, 


y) 


1 Xj-Vj 

Att |x — yP ’ 


( 4 ) 

( 5 ) 
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where 5ij is the Kronecker delta. More precisely, the velocity vector u(x) and pressure p(x) are 
given by 

3 3 

Ui{^) = p(x) = ^Pj(x,y)/j-. (6) 

i=i i=i 

The stresslet (or double layer kernel) for the Stokes equations describes the velocity induced by 
an infinitesimal displacement g = {gi, 92 , 93 ) (sometimes called the double force source strength) at 
a point y with orientation vector 1 / = [vi, 122 , 1 ^ 3 )■ It is given by: 


Tijki^, y) 


3 {xj - yi){xj - yj){xk - Vk) 

dvr |x — yp 


hj,k = 1,2,3, 


( 7 ) 


njfc(x,y) = ^ 


6jk Hxj - yj){xk - Vk) 


+ 


|x-y|'" |x-y|'" 

and the corresponding formulas for the velocity and pressure at an arbitrary point x are 


( 8 ) 


Tijk{y^,y)vkgj, p{^) = EE njfc(x,y)zyfc5'j. 


(9) 


j=i k=i 


j=i k=i 


2.1 Stokes flow in a half-space, Blake’s formula 

Suppose now that a force vector f is applied to a viscous fluid in the upper half-space (xa > 0). 
Then the corresponding Stokeslet-induced velocity field fails to satisfy the no-slip condition Q. 
In order to annihilate the velocity field while satisfying the homogeneous Stokes equations in the 
upper half-space, [1] proposed the following image structure: 

‘S'*^(x,y) = 5ij(x,y)-5ij(x,y^)-b2y|5,^(x,y^)-2y35,^^(x,y^), (10) 

where y^ = {yi, 92 , — 93 ) is the reflected image location. Here, is a modified source doublet given 
by 

= (11) 

and is a modified Stokes doublet, given by 

S.f(x,y) = (i-2i,-3)^^i^)T, (12) 

Similar, but more involved, decompositions for the Stokes doublet and stresslet in a half-space are 

given in and [20j. Note that the computation of the modified source doublet S[j requires the 

evaluation of three distinct harmonic dipole fields. 


3 A new image formula 

In this section, we derive a simpler image structure, using the Papkovich-Neuber representation 
which involves only a single harmonic function, yfi as above, will denote the image location 
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(yi)y 2 ) —ys)- Now, however, we define reflected single force, double force, and double force orien¬ 
tation vectors by negating their third components: 

= (/i,/2,-/3), = ( 51 , 52 ,- 53 ), = { 1 ^ 1 ,1^2,- 1 ^ 3 )- (13) 

We also make use of the harmonic potential due to a unit strength charge. 


G^(x,y) = 


1 


47r|x — y| ’ 

the harmonic potential due to a unit strength dipole with orientation vector i', 


d 


1 


i=l 


dyi 47r|x - y| 


(14) 


(15) 


and the harmonic potential due to a unit strength quadrupole with orientation vectors u and k. 


G'^ [u, k] (x, y) = 

i=l j=l 


92 


1 


dyidyj 47r|x - y| 


(16) 


Finally, recall that, through the Papkovich-Neuber representation Q, the harmonic function (/)(x) 
induces velocity and pressure fields that can be written in component form as 


0 d 

Mi(x) = aia—(/)(x) - (ii 3 </>(x), p(x) = 2-—(/>(x) 

C/Xi 


(17) 


Of particular note is the fact that at X 3 = 0, the only non-zero velocity component is 


«3(x) = -(/>(x). 


(18) 


This suggests a simple two-step strategy. First, annihilate the tangential components of the 
velocity field induced by a Stokeslet (single layer) or stresslet (double layer) kernel. It is easy 
to see that this can be accomplished by subtracting the influence of a reflected single force = 
(/i,/ 2 ,-/ 3 ), or double force = ( 51 , 52 ,- 53 ),!''^ = (z^i, i^ 2 ,-z^s) located at the image point y^, 
respectively. It remains only to match the remaining non-zero normal component U 3 , which we will 
do by a judicious choice of the Papkovich-Neuber potential </>(x). 


3.1 The Stokeslet correction 

Following the discussion above, let us write out in more detail the velocity field v remaining after 
subtracting the image Stokeslet at from the original Stokeslet at y (Fig. 1); 

3 3 


j=i j=i 


(19) 


At the interface X 3 = 0, a straightforward computation yields 

fi 


Ul(x) = 0, U 2 (x) = 0, U 3 (x) = - 


1 


dvr |x — 


ys Xj yj j 

lx — y^p ^ 
j=i I •>' I 


( 20 ) 


4 








Figure 1: The original source location is in the upper half-space at (yi, 2 / 2 ) 2 / 3 ) and the reflected 
image source is at ( 2/1 j2/2)— 2 / 3 )• For target points (xi,X 2 , 0 ) that he on the half-space boundary, 
the distances r and R from the source and image are the same, simplifying the computation of v 
in pOl). 


By inspection of (20), it is clear that the harmonic potential (j) required to cancel the non-zero 
normal component U 3 at the interface is that induced by a charge of strength /g and a dipole with 
orientation vector of strength 2/3 located at y^. Thus, the velocity field satisfying the desired 
no-slip boundary condition ([^ can be expressed as: 

uV (x) = uf{x) - uf(x) - uf (x), (21) 

where and are the velocity fields induced by the original free-space Stokeslet and the reflected 
image force vector, respectively: 


u. 


= '^Sij{x,y)f^ 
i=i 




Here, uf is the Papkovich-Neuber correction 


d 


u. 


B 




( 22 ) 


Ui (x) = X3 — (j){x) - 6i3(j){x), 


(23) 


where the harmonic potential ([> is that due to a simple charge and dipole, both located at y^: 

4>{x) = fi G^(x,y^) + y 3 G^[f'](x,y^). (24) 


3.2 The stresslet correction 

We state the image structure for the stresslet (the double layer kernel) in the form of a theorem. 
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Theorem 1. Let a double force g with orientation vector v be located at y, resulting in the free-space 
velocity field 

3 3 

Tijk{x,y)ukgj. (25) 

j=i k=i 

Then the corresponding velocity satisfying the no-slip boundary condition is given by 
3 3 3 3 

y)^k9j EE Tijk{x, 

j=l k=l j=l k=l 


t/Wkd] - 


X3^4>^{x) - discfy^ix) 

UXi 


(26) 


where 

(j)^{x) = 2{u^ • g^)G^[h]{x, i/) + 2y3G*3[i/^,g^](a;, t/), 

with h = (0,0,1). 

Proof. Note first that the symmetric part of the Stokes doublet is given by 


^l(x,y) = ^ 


dvr 


{xi-yi) 3{xi - yi){xj - yj){xk - yk) 

-+- 


- yf 


|x - y| 


Furthermore, the tangential components of the velocity held 

3 3 3 3 

Vi{^) = 

j=l k=l j=l k=l 


(27) 


(28) 


(29) 


are zero when X 3 = 0. Thus, the harmonic potential annihilating at the interface X 3 = 0 is 
that due to a single quadrupole source located at with orientation vectors and strength 

22 / 3 : 

(()^(x) = 2y3G«[i.^g^](x,y^). (30) 


It is straightforward to check that the stresslet is simply related to the symmetric part of the Stokes 
doublet 


7iifc(x,y) 


3 {Xj - yi){xj - yj){xk - yk) 
dvr |x — yp 


1 Xi-yi 
A I I 3 ^ ^ 

47r |x — 


+ ^ifc(x,y)- 


(31) 


Thus, it remains only to annihilate the velocity held induced by the hrst term on the right-hand 
side of (31) and its rehected image, given by 


fpx) 


1 

dTT 


3 3 


EE 


Vi 

X - yP 


djkk'kdj 


1 

d-TT 


3 3 

EE 

j=l k=l 


Xi - Vi 


^jk’^igj ■ 


X - y 


/|3 


(32) 


At X 3 = 0, the tangential components vanish and the normal component is easily computed to be 
2 ( 1 /^ • g^)G^[h](x, y-^), where h = (0,0,1). The desired result follows. □ 


6 













3.3 The rotlet correction 

Similar representations are easily derived for other fundamental solutions, such as the rotlet - the 
antisymmetric part of the Stokes doublet; 


^ffc(x,y) = iz 


1 


'-ijk 


dvr 




3 


|x-y|-^ |x-y| 

The corresponding velocity satisfying the no-slip boundary condition Q is given by 
3 3 3 3 


w = EE 


d 


j=l k=l j=l k=l 

where the Papkovich-Neuber correction (f)^ is due to two dipoles: 

= -2i^|G^[g^](x,y^) -h 25 |G^[i/^](x,y^). 




(33) 


(34) 


(35) 


3.4 The Stokes doublet correction 

Finally, the Stokes doublet is the sum of its symmetric and antisymmetric parts 


Tlfki^,y) = iz 


1 


dvr 


(xi - Vi) 


|x - y 


3 ^:ik 


^ + i^k-vk) ^ 


|x - yi 


|x - y 


+ 


3{xi - yi){xj - yj){xk - yk) 


|x - y| 


3 


= 7i‘vfc(x,y) + 7)ffc(x,y). 


By combining (30) and (35), we obtain the image formula for the Stokes doublet 

3 3 


3 3 


D,W 

Ui 


w = EE TFjk{^,Y)yk9j - J2Y.^Sk(^^y^yk9: 

j=l k=l j=l k=l 

where the Papkovich-Neuber potential is given by 

</.^(x) = -2v.|G^[g^](x,y^) + 25 |G^[*.^](x,y^) + 2v/3GQ[^.^g^](x,y^). 


(36) 


(37) 


(38) 


4 Conclusions 

We have derived very simple image formulas for Stokes flow in a half-space induced by any of the 
standard fundamental solutions - the Stokeslet, stresslet, rotlet, and Stokes doublet. In each case, 
all that is required is a reflected fundamental solution and a Papkovich-Neuber correction based on 
a single harmonic potential. 

The velocity (and pressure) due to the “direct” and reflected fundamental solutions can be 
computed together with any software that handles Stokeslets, stresslets, etc. in free space. Fur¬ 
thermore, the Papkovich-Neuber potential requires only the evaluation of a single additional har¬ 
monic function — itself requiring only software for free space harmonic sources, dipoles and 
quadrupoles. Many efficient schemes exist for these various steps, such as those described in 

linillllEElEZlIlEaiislEelEllESlEllIS]. 
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A Extension to the two-dimensional problems 


The two-dimensional Stokes flow representations in a half-space can be derived similarly. They 
lead to identical Papkovich-Neuber corrections with corresponding charge, dipole, and quadrupole 
potentials replaced by their two-dimensional equivalents. Similar to the three-dimensional case, 
for any x G IR^, the velocity field satisfying the no-slip boundary condition in a half-plane is 
composed of three terms: 

u’^(x) = u^(x)-u^(x)-u^(x), (39) 

where the first term is the velocity field induced by the free-space Green’s function, the second 
term is the reflected image about the plane wall annihilating the tangential velocity component, 
and the third term is a Papkovich-Neuber correction term in the following form: 


u‘^(x) = X 2 V 4 >{x) 


0 


p(x) = 2 


dx2 


(40) 


The correction potentials for various fundamental solutions are 

<))(x) = fi G^(x,y^) + y 2 G^[f'](x,yO, (41) 

for the two-dimensional Stokeslet, and 

(/)'^(x) = 2(r/^ . g^)G^[h](x,y^) + 2y2G«[i/^ g'](x, y^), (42) 

0^(x) = -2i/^G^[g^](x,y^) -F2g|G^[iv^](x,y^), (43) 

</.^(x) = -2z.|G^[g^](x,y^) + 25|G^[*.^](x,y^) + 22/2G'Q[^.^g^](x,y^), (44) 

for the two-dimensional stresslet, rotlet, and Stokes doublet, respectively, where the images of the 
source, single force, double force, and double force orientation vectors are 

= {yi, -y 2 ), = (/i, -/ 2 ), = (gi, - 92 ), = ( 1 ^ 1 , - 2 ^ 2 ), (45) 


with the orientation vector h = (0,1), and the free-space Laplace Green’s functions in two dimen¬ 
sions are given by 




1 ^ 8 
(x,y) = log|x-y|, G^[i/](x,y) = y), 


2 2 


2=1 


and G^[i/,K](x,y) = y) 


i=i j=i 


(46) 


(47) 


B Extension to linear elasticity kernels 

The single layer kernel for linear isotropic elasticity in is given by Kelvin’s solution HZ! 

1 


Uij{x,y) = 


S-Kfl 




- y| 


|x - y| 


, bj = 1,2,3, 


(48) 











where A, fi are Lame’s parameters and a = (A + /i)/(A + 2;u). It is easy to see that the displacement 
field 

3 3 

j=i i=i 


satisfies the no-displacement boundary condition u(x) = 0 at X 3 = 0, where the Papkovich-Neuber 
correction is 


d 

up(x) = ax 3 —(f{x) - (2 - a)5i34>{x). 

(50) 

<^(x) = ^G^(x,y^) + ^G^[f^](x,y^). 

/i p, 

(51) 
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